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Abstract 

We report the observation of two narrow structures at 10610 MeV/c^ and 10650 MeV/c^ in the 
7r^T(n>S') (n = 1,2,3) and n^hb{mP) (m = 1,2) mass spectra that are produced in association 
with a single charged pion in T(55') decays. The measured masses and widths of the two structures 
averaged over the five final states are Mi = 10608.4 ± 2.0MeV/c^ Ti = 15.6 ± 2.5 MeV and M2 = 
10653.2 ± 1.5 MeV/c^, T2 = 14.4± 3.2MeV. Analysis favors quantum numbers of /^(J^)=l+(1+) 
for both states. The results are obtained with a 121.4 fb~^ data sample collected with the Belle 
detector near the T(56') resonance at the KEKB asymmetric-energy e^e~ collider. 

PACS numbers: 14.40.Pq, 13.25.Gv, 12.39.Pn 
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INTRODUCTION 



Analyses of decays of the T(55') resonance to non-BsBs final states have produced several 
surprises. Recently the Belle Collaboration reported observation of anomalously high rates 
for T{5S) T{nS)7r+7r- {n = 1,2,3) [l] and T(5S) ^ h{mP)7r+7T- (m = 1,2) 
transitions. If the T{nS) signals are attributed entirely to T(5S') decays, the measured 
partial decay widths r[T(5S') — ?> T(r2S')7r+7r^] ~ 0.5 MeV are about two orders of magnitude 
larger than typical widths for dipion transitions among the four lower T{nS) states. The 
rates for the T(55') — ?■ h},{'mP)'K^ -k^ transitions are found to be comparable with those for 
T(5S') — )■ T(riS')7r+7r^, indicating that processes requiring a heavy-quark spin flip - i.e., in 
hh{mP) production - do not seem to be suppressed. These unexpected observations indicate 
that exotic mechanisms could be contributing to the T{5S) decays. 

We report preliminary results of detailed studies of the three-body T{5S) — )■ T(nS')7r+7r~ 
(n = 1,2,3) and T(55') — )■ hb{mP)Tr~^n~ {m = 1,2) decays. Results are obtained with a 
121.4 fb~^ data sample collected near the peak of the T(5S') resonance (y^ ~ 10.865GeV) 
with the Belle detector at the KEKB asymmetric energy e~^e~ collider j3|. 



BELLE DETECTOR 



The Belle detector is a large-solid-angle magnetic spectrometer that consists of a 4-layer 
silicon vertex detector (SVD), a 50-layer central drift chamber (CDC), an array of aerogel 
threshold Cherenkov counters (ACC), a barrel-like arrangement of time-of- flight scintillation 
counters (TOF), and an electromagnetic calorimeter (ECL) comprising CsI(Tl) crystals 
located inside a superconducting solenoid coil that provides a 1.5T magnetic field. An iron 
fiux-return located outside of the coil is instrumented to detect Kr mesons and to identify 

n 

muons (KLM). The detector is described in detail elsewhere 

For charged hadron identification, the dE/dx measurement in the CDC and the response 
of the ACC and TOF are combined to form a single likelihood ratio. Electron identification is 
based on a combination of dE/dx measurement, the response of the ACC, and the position, 
shape and total energy deposition of the shower detected in the ECL. Muons are identified 
by their range and transverse scattering in the KLM. 

We use a GEANT-based Monte Carlo (MC) simulation |5[ to model the response of the 
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detector and determine the acceptance. The MC simulation includes run-dependent detector 
performance variations and background conditions. 

ANALYSIS OF T(55) T(15, 25, 35)7r+7r- 

To select T(5S') — )■ T(nS')7r+7r^ candidate events we require the presence of a pair of muon 
candidates with an invariant mass in the range of 8.0 GeV/c^ < M(/i "•"/!") < 11.0 GeV/c^ 
and two pion candidates of opposite charge. These tracks are required to be consistent with 
coming from the interaction point: we apply requirements dr < 0.5 cm and \dz\ < 2.5 cm, 
where dr and \dz\ are impact parameters perpendicular to and along the beam axis with 
respect to the interaction point, respectively. We also require that none of the four tracks 
is positively identified as an electron. No additional requirements are applied at this stage. 

Candidate T{5S) — )■ T(nS')7r+7r~ events are identified by the invariant mass of the fi~^^~ 
combination and the missing mass MM(7r+7r~) associated with the vr+vr" system calculated 
as 

MM(7r+7r-) = ^ (E,^^^ - E^.^^y ~ p*X., (1) 

where Ec.m. is the center-of-mass (cm.) energy and E*+^_ and are the energy and 

momentum of the 7t~^7t~ system measured in the cm. frame. The two-dimensional distribu- 
tion of M{fi'^fi~) versus MM(7r+7r~) for all the preselected candidates is shown in Fig. [TJ 
Events originated from T(5S') decay fall within the narrow diagonal band (T(5S') signal 
region) defined as |MM(7r+7r") -M(/i+yU")| < 0.2 GeV/c^ (see Fig. [1]). Clustering of events 
around nominal T{nS) mass values are clearly visible on the plot. 

The amplitude analyses of the three-body T(55') — )■ T(n5')7r+7r^ decays that are reported 
here are performed by means of unbinned maximum likelihood fits to two-dimensional Dalitz 
distributions. 

Before fitting the Dalitz plot for events in the signal region, we determine the distri- 
bution of background events over the Dalitz plot using events in the T{nS) mass side- 
bands that are refitted to the nominal mass of the corresponding T{nS) state to match 
the phase space boundaries. Definitions of the mass sidebands are given in Table [T] and 
the corresponding sideband Dalitz plots are shown in Fig. [2l where M{T{nS)Tc)^in and 
M(T(n^)7r)max are defined as M{T{nS)7T)^i^ = mm{M{T{nS)7T+), M{T{nS)n-)) and 
M(T(n5)7r),„ax = max(M(T(n5)7r+), M(T(n5)7r-)), respectively. 
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FIG. 1: A scatter plot of all of the T(5S') — t- T(nS')7r+7r candidate events that pass the initial 
selection criteria. The area between the two diagonal lines is the T{5S) signal region. 



CI 

a 
- 110 



t 106 
104 



/ ■ ■ 


(a) 








(■■■■.." 




-y- ■■ ■■■ :.■ ; 
-f ■ . ' 




■ 1 '• ' ■ ■ ■ . 
. - ■■ ' 

■ , , , , 1 , , , 


. , , . 1 , , . , 1 



il12 



(b) 



(1) 
o 



0.2 0.4 0.6 

(7t*7c") , GeV^/c" 



I- 



(c) 



(7U*7C ) , Ge\^/c 



FIG. 2: Dalitz plots for T(n5)^+^- events in the (a) T{1S); (b) T(25); (c) T(35) sidebands. 
It is evident in the sideband Dalitz distributions, shown in Fig. |2l that there is a strong 



TABLE I: Definition of the sidebands for T(nS')7r+7r signals. The sidebands are defined in terms 
of regions in the spectrum of the missing mass associated with the 7r"^7r~ system. 



Final state 


T{1S)tt+tt- 


T{2S)tt+tt- 


T{3S)tt+tt- 


Sideband 
region, GeV/c^ 


9.38 < MM(7r7r) < 9.43 
9.48 < MM(7r7r) < 9.53 


9.95 < MM{-ktt) < 10.00 
10.05 < MMi'KTT) < 10.10 


10.28 < MMi'KTT) < 10.33 
10.38 < MMi'KTT) < 10.43 
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concentration of background events in the very low ti^ti~ invariant mass region; these are 
due to photon conversion in the innermost parts of the Belle detector. Because of their 
low energy, these conversion electrons are poorly identified and pass the electron veto re- 
quirement. We exclude this high background region by applying the requirements on the 
TT+TT^ invariant mass given in Table HTl For the remainder of the Dalitz plot the distribution 
of background events is assumed to be uniform. The variation of reconstruction efficiency 
across the Dalitz plot is determined from MC simulation. The fraction of signal events in 
the signal region for each of the three T(?2S')7r+7r^ final states is determined from a fit to 
the corresponding MM(7r+7r^) spectrum using a Crystal Ball function [7] for the T signal 
and a linear function for the combinatorial background component. Results of the fits are 
shown in Fig. [3] and are summarized in Table [ITl 

Figure m shows Dalitz plots of the events in the signal regions for the three decay channels 
under study. In all cases, two horizontal bands are evident in the T{nS)Ti system near 
10.61 GeV/c^ (~ 112.6 GeVVc^) and 10.65 GeV/c^ (~ 113.3 GeVVc^). In the following we 
refer to these structures as Zf,(10610) and Zfc(10650), respectively. 

We use the following parameterization for the T(5S') — )■ T(nS')7r+7r~ three-body decay 
amplitude: 

M(si, S2) = Ai{si, S2) + A2{si, S2) + Af„ + Af^ + Anr, (2) 

TABLE II: Results of the fits to the MM{iT^7r~) missing mass distributions. 

Final state T(15)7r+7r- T{2S)tt+tt- T{3S)tt+tt- 

M(7r+7r-) signal 0.20 < M(7r+7r") 0.16 < M(7r+7r") 0.10 < M(7r+7r") 
region, GeV/c^ 

MM{tt+7t-) signal 9.43 < MM{tt+tt-) 10.00 < MM{tt+tt-) 10.33 < MM(7r+7r-) 

region, GeV/c^ MM{tt+tt') < 9.48 MM{tt+tt-) < 10.05 MM{tt+tt-) < 10.38 

Peak position, 9455.3 ± 1.1 10019.1 ± 1.7 10350.2 ± 2.3 
GeV/c2 

Resolution, 7.31 6.74 6.33 

MeV/c2 

Number of events 1581 1926 689 

Signal fraction 0.902 ± 0.015 0.937 ± 0.007 0.917 ± 0.010 
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FIG. 3: Distribution of missing mass associated with the 7r"'"7r~ combination for T(5«S') — >■ 
T(nS')7r+7r^ candidate events in the (a) T(IS'); (b) T{2S); (c) T{3S) mass regions. Vertical 
Hnes define the corresponding signal regions. 
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FIG. 4: Dalitz plots for T(n5)^+7r" events in the (a) T{1S); (b) T(25); (c) T(35) signal regions. 
Dalitz plot regions to the right of the vertical lines are included in the amplitude analysis. 

where Si = M'^{T{nS)Tc^), S2 = M^(T(n5')7r^). Here we assume that the dominant con- 
tributions come from amplitudes that conserve the orientation of the spin of the heavy 
quarkonium state and, thus, both pions in the cascade decay T{5S) — ?■ Z^tt — t- T(r2S')7r+7r~ 
are emitted in an S-wave with respect to the heavy quarkonium system. As will be shown 
later, angular analyses support this assumption. Consequently, we parameterize the ob- 
served Zb(10610) and Z;,(10650) peaks with an S-wave Breit-Wigner function 

r 



BW{s,M,T) 



(3) 



M2 - s - iMT ' 

where we do not consider possible s— dependence of the resonance width F. To account 
for the possibility of T(55') decay to both Z^tt~ and Z~tt'^, the amplitudes Ai and A2 
are symmetrized with respect to tt'^ and vr" transposition. Taking into account isospin 
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symmetry, the resulting amplitude is written as 

Akisi, S2) = aue^'^ {BW{s^, M^, r^) + BW{s2, M^, r^)), (4) 

where the masses and the widths {k = 1, 2) are free parameters of the fit. Due to the 
very limited phase space available in the T{5S) T(3S')7r+7r~ decay, there is a significant 
overlap between the two processes T{5S) — t- Z^'k^ and T(5S') — t- Z^t^^ ■ We also include Af^ 
and Aj^ amplitudes to account for possible contributions in the 7r"*'7r~ channel from /o(980) 
scalar and /2(1270) tensor states, respectively. Inclusion of the /o(980) state is necessary 
in order to describe the prominent structure in the M(7r+7r~) spectrum for the T(lS')7r+7r~ 
final state around M(7r"''7r~) = 1.0 GeV/c^ (see Fig. [5]). We also find that the addition of the 
/2(1270) gives a better description of the data at M{t\^i\^) > l.OGeV/c^ and drastically 
improves the fit likelihood values. We use a Breit-Wigner function to parameterize the 
/2(1270) and a coupled- channel Breit-Wigner (Flatte) function ^] for the /o(980). The 
mass and the width of the /2(1270) state are fixed at their world average values [61]; the 
mass and the coupling constants of the /o(980) state are fixed at values defined from the 
analysis of 5+ ^ K+n+n-: M(/o(980)) = 950 MeV/c^, g^^ = 0.23, qkk = 0.73 

Following suggestions given in Refs.|lo| and references therein, the non-resonant ampli- 
tude Anr has been parameterized as 

ANR = are^^?'' + a-e*^"''s3, (5) 

where S3 = M^(7r^7r^) (S3 is not an independent variable and can be expressed via Si and 
S2 but we use it here for clarity), , 02^, ^i"^ and are free parameters of the fit (with an 
exceptions of the T(3S')7r+7r^ channel as described below). 
The logarithmic likelihood function C is then constructed as 

C = -2J2 log(/sig5(si, S2) + (1 - /sig)i?(si, S2)), (6) 

where S'(si,S2) = |M(si,S2)P folded with the detector resolution function (2 MeV/c^ for 
M{T{nS)TT^) when both T{nS) and T{5S) mass constraints are used; the M(7r+7r^) resolu- 
tion is not taken into account since no narrow resonances are observed in the tt+tt^ system), 
-B('Si,S2) = 1 and fsig is a fraction of signal events in the data sample. Both S'(si,S2) and 
B{si, S2) are corrected for reconstruction efficiency. 

In the fit to the T(lS')7r+7r^ sample, the amplitudes and phases of all of the components 
are allowed to fioat. However, in the T(25')7r+7r^ and T(35')7r+7r^ samples the available 
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FIG. 5: Comparison of fit results (open histogram) with experimental data (points with error bars) 
for events in the T{1S) (top row), T{2S) (middle row), and T{3S) (bottom top) signal regions. 
The hatched histogram show the background component. 

phase space is significantly smaller and contributions from the /o(980) and /2(1270) are not 
well defined. Thus, in the fit to the T(2S')7r"^7r~ and T(3S')7r^7r~ signal samples, we fix the 
amplitudes and relative phases of these components to the values measured in the fit to the 
T(lS')7r^7r^ sample. Moreover, in the fit to the T(3S')7r^7r~ sample, we also fix the and 62'^ 
parameters of the A^^ amplitude. Possible effects of these assumptions are considered while 
determining the model-dependent uncertainty. Results of the fits to T(5S') — )■ T(nS')7r+7r~ 
signal events are shown in Fig. [5l where one-dimensional projections of the data and fits 
are compared. To combine and events we plot T{nS)7c mass distributions in terms 
of M(T(nS')7r)min and M(T(nS')7r)max; fits are performed in terms of M(T{nS)TT'^) and 
M{T{nS)TT^). Results of the fits are summarized in Table IVll We try various alternative 
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models to parameterize the decay amplitude as described in the systematic uncertainty 
section. The combined statistical significance of the two peaks exceeds 10 a for all tested 
models and for all T(nS')7r+7r^ channels. 

The main sources of systematic uncertainties are found to be: 

• Dependence of reconstruction efficiency on the Dalitz plot variables has been deter- 
mined from MC simulation of signal events uniformly distributed over the Dalitz plot. 
Since the pions from the T(55') — ?■ T(r2S')7r+7r^ decays have low momenta, the effi- 
ciency determined from MC might be systematically shifted by up to 5-10%. The 
effect on the parameters of the observed peaks is estimated by varying the efficiency 
at the edges of the Dalitz plot by ±10% and refitting the data. 

• To estimate the effect of the parameterization used to fit the data on the parameters 
of the Zb peaks we fit the data with modifications of the nominal model described in 
Eq. |2l In particular, we vary the M(7r+7r^) dependence of the non-resonant amplitude 
Anr, try to include a D-wave component into Anr, etc. The variations in the extracted 
Zh parameters determined from fits with modified models are taken as estimates of 
the model uncertainties. 

• The uncertainty in the cm. energy leads to an uncertainty in Dalitz plot boundaries. 
This effect is particularly important for the T(3S')7r+7r^ channel. To estimate the 
associated effect on the Zb parameters, we generate normalization phase space MC 
samples that correspond to E^m. ± 3 MeV, where Ecm. is the nominal cm. energy. 

Systematic effects associated with uncertainties in the description of the combinatorial 
background is found to be negligible. The results of our studies on systematic uncertainties 
are summarized in Table IIIII 

ANALYSIS OF T(5S') hb{lP,2P)7T+7T- 

hb{lP) and hb{2P) have been recently observed by Belle in the decay T{5S) — )■ 
hb{'mP)7r~^7i~ [21]. Here we study resonant structure of these decays. Because of high back- 
ground Dalitz plot analysis is impossible with current statistics, therefore we study the one 
dimensional distributions in M{hb{mP)Tx). We use the same selection requirements and 
reconstruction procedure as in Ref. j2|. 
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TABLE III: Summary of dominant sources of systematic uncertainties. Uncertainties for masses 
and widths are given in MeV/c^ and in degrees for the relative phase. Quoted numbers are for 
T(1S')/T(2S')/T(3S') channels, respectively. 



Parameter 


Efficiency 


ivioctei 




lotai 


M(Z6(10610)) 

r(Zb(io6io)) 

M(Zb(10650)) 


±1/±1/±1 
±1/±1/±1 
±1/±1/±1 


I 9/+0.3 /+2.4 
='=^/-1.4/ 1.2 


±l/±l/±2 
±l/±2/±3 
±1/±1/±1 


±2/l^/+5 

±2/tI/±4 
±2/±2/±2 


r(Zfe(10650)) 


±1/±1/±1 




±l/±3/±l 




Rel. phase 


±3/t^/±5 


+0 /+11 /+21 
-50/ -6 /-58 


±3/±6/±9 


+5 /+14/+23 
-50/ -9 /-59 


Rel. amp. 


±0.01/ ±0.02/ ±0.02 


+0.08 /+0.03 /+0.14 
-0.0 /-O.O /-0.02 


±0.02/ ±0.02/ ±0.04 : 


+0.09 /+0.04 /+0.15 
-0.03/ -0.03/ -0.0.05 



The T(55') — ?■ hh{mP)7T~^7[~ decays are reconstructed inclusively using the missing mass 
of the Tc'^TT^ pair. We select 7r+ and tt~ candidates that originate from the vicinity of the 
interaction point {dr < 3 mm, \dz\ < 2 cm) and are positively identified as pions. We 
reject tracks that are identified as electrons. The continuum e~^e~ qq {q = u, d, s, c) 
background is suppressed by a requirement on the ratio of the second to zeroth Fox- Wolfram 
moments R2 < 0.3 fll)]. 

We define the M{hb{mP)'n'~^) as a missing mass of the opposite sign pion, MM(7r^). We 
measure the yield of signal decays as a function of the MM{tt^) by fitting the MM(7r+7r^) 
spectra in the bins of MM{7c^). We combine the MM^n^n^) spectra for the corresponding 
MM{7i~^) and MM(7r~) bins and we use half of the phase-space to avoid double counting, 
as shown in Fig. [61 

In the studies of the T{5S) — > /i;,(lP)7r+7r^ channel we subdivide the interval 
10.4 MeV/c^ < MM{tt) < 10.72 MeV/c^ into 32 bins and perform a fit to the MM{n+TT-) 
spectrum for each MM{t[) bin. The fit function consists of four components: the /ib(lP) 
signal, the T{2S) signal, a reflection from the T(3S') — t- T(lS')7r+7r^ decay, and combinato- 
rial background. The shapes of the signal components and T(35') reflection is determined 
from the analysis of exclusive /i^/i^Tr+vr^ data as described in Ref. [2I. The signals are 
parameterized by a Crystal Ball function to accommodate tails due to initial state radiation 
of soft photons that accounts for about 8% of the signal yield. The resolution width of the 
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FIG. 6: The MM{tt~) versus MM{tt+) distribution for the T(55) hb{lP)7r+TT- (left) and 
T(5S') — )• /ife(2P)7r+7r^ (middle) events simulated according to the phase-space MC model. Lines 
indicate the values of cuts used to avoid double counting after symmetrization relative to vr"'" and 
7r~ . Right: the MM(tt) versus MM(7r+7r^) distribution for inclusive K'^ — )■ ir^ir" signal in the 
fast MC simulation. The red rectangle shows the fit region for the hb{2P)'K^ii^ analysis. 

T(25) {K{IP)) is a = 6.5MeV/c2 (G.SMeV/c^). The T(35) ^ T(15)7r+7r- reflection is 
described by a single Gaussian function with the width of a = 18MeV. The ratio of the 
T(3S') — )■ T(lS')7r+7r^ and T(2S') yields is measured from /i+yU^vr+Tr^ data and corrected 
for the branching fractions for T(15') — )■ /i"'"yU~ and T{2S) — t- yU^/x^ decays. In the fit to 
the inclusive MM(7r+7r^) spectra, this ratio is allowed to float within the uncertainties of 
the exclusive fi~^ fi^TT^n^ measurements. The combinatorial background is parameterized by 
a Chebyshev polynomial. The polynomial order decreases monotonically from 10-th order 
for the first bin to 6-th order for the 11-th bin. The background shape in the first MM{tt) 
bins is complicated due to the proximity of the kinematic boundary, therefore a high-order 
polynomial is used. 

Results of fits for the h},{lP) and T(25') yield as a function of MM{tt) are shown in Fig. [71 
where the T(2S') yield measured from the analysis of the exclusive T(5S') yU+yU^7r+7r~ 
data and renormalized to the total number of events in inclusive spectrum is also shown for 
comparison. The T(2S') yields for inclusive and exclusive measurements agree well. The 
hi,{lP) yield shown in Fig. [7] exhibits a clear two-peak structure without any significant 
non-resonant contribution. The peak positions are consistent with those for the ^{,(10610) 
and Zf,(10650) observed in the T{nS)iT^ final states. 

We perform a fit to the MM(7r) distribution. We assume that spin-parity for both 
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FIG. 7: Left: the yield of T{2S) as a function of the MM{'n) measured using the inclusive data 
(points with error bars) and exclusive /x^/U^vr^Tr^ data (histogram). Middle: the yield of the 
h{,(lP) as a function of MM{tt) (points with error bars) and results of the fit (histogram). Right: 
the yield of the hb{2P) as a function of MM{'k) (points with error bars) and results of the fit 
(histogram). 

Z{,(10610) and ^{,(10650) is = therefore in the fit function we use a coherent sum of 
two P-wave Breit-Wigner amplitudes; we add also a non-resonant contribution. 



f = A \BW{s, Ml, Ti) + ae''^BW{s, M2, T2) + be'^^l^ 



(7) 



Here ^/s = MM{it); the variables A, M^, [k = 1,2), a, 0, b and are floating in the fit; 
^ is a phase-space factor, p (g) is the momentum of the pion originating from the T{5S) 
(Zb) decay measured in the rest frame of the corresponding mother particle. 
The P-wave Breit-Wigner amplitude is expressed as 

q/^sF 



BW{s,M,T) 



M2 



s-iMT(s) 



(8) 



Here F is the P-wave Blatt-Weisskopf form-factor F = ^ |12|. go is daughter mo- 

mentum calculated assuming pole mass of its mother, R = 1.6GeV~^; F(s) is the energy- 
dependent width, F(s) = F(^)^-^p2. The function / is convolved with the detector resolu- 
tion function, is integrated over the lOMeV/c^ wide bin and is corrected for reconstruction 
efficiency. The detector resolution is parameterized by a single Gaussian function with 
a = 5.2 MeV/c^ as determined from MC simulation. The result of the fit is shown in Fig. [7] 
and is summarized in Table IVII The non-resonant amplitude is found to be consistent with 
zero, b = 0.03 ± 0.04. The confidence level of the fit is 81%. We find that the hypothesis of 
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two resonances is favored over the hypothesis of a single resonance (no resonances) at the 
7.4 (T (17.9 (t) level. 

The systematic uncertainty is estimated by varying (within ±2) the order of the Cheby- 
shev polynomial and repeating the fit to the data. We also perform fits with the constraint 
on the normalization of the T(3S') — t- T(lS')7r+7r^ reflection released. We shift the MM{it) 
binning by half the bin-size. To estimate the uncertainty due to fit model, we try to fit 
the MM(7r) spectrum with the non-resonant component fixed at zero. We vary the Blatt- 
Weisskopf parameter R from to 5GeV~^. As found in Ref. [2I, the detector resolution 
as determined from a MC simulation could be underestimated by as much as 5-10%. To 
account for this effect we repeat the fit to the MM (vr) spectrum with a corrected resolu- 
tion function. The maximal change of signal parameters for each source is considered as a 
systematic uncertainty. We find that the systematic uncertainties associated with releasing 
the constraint on the T{3S) — ^ T(lS')7r+7r~ reflection and varying R are negligible. Finally, 
we find about IMeV/c^ deviations of the T(IS'), T{2S) and T(3S') peak positions relative 
to the PDG values when the T{nS) states are reconstructed inclusively j^. These small 
deviations could be due to the local variations of the background shape that are not fully 
described by the polynomial. Consequently, we add an additional uncertainty of IMeV/c^ 
to all mass measurements. All contributions are added in quadrature to obtain the to- 
tal systematic uncertainty. The summary of systematic uncertainties study is presented in 
Table IIV[ The minimal level at which the two-resonance hypothesis is favored over the 
one-resonance (no-resonance) hypothesis for all considered variations is 6.60" (16.0 a). 

The analysis of the hb{2P) tt+tt^ final state follows the same strategy as that for the 
hf,{lP) n^n^ . In this case we require MM(7r) > 10.57 GeV/c^ to avoid double counting of 
signal events. We also require MM(7r+7r-) < 10.34 GeV/c^ and MM{7i) < 10.7GeV/c2 to 
minimize influence of the reflections from events where the two charged pions originate from 
a ^ TT+TT- decay (see Fig. [6]). We subdivide the interval 10.57 GeV/c^ < MM{n) < 
10.7GeV/c^ into 13 bins and perform fits to the MM(7r+7r^) spectra for each MM{tt) 
bin. The flt function consists of three components: the hb{2P) signal, reflection from the 
T{2S) — 7- T(lS')7r+7r^ decay and combinatorial background. The T(3S') signal is outside 
the mass range. The shapes of the signal and reflection are determined from the analysis of 
the exclusive T(55') — ?■ fi~^ fi^n^TT^ data. The yield of the reflection for each MM(7r) bin is 
determined from the exclusive /i+yU^Tr+vr^ data and is normalized to the ratio of inclusive 
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TABLE IV: Systematic uncertainties in the mass and width of the ^{,(10610) (index 1) 
and 10650) (index 2), their relative normalization factor a and phase <j) for the 
/i;,(lP)7r+7r^//i;,(2P)7r+7r~ decay modes. 





Chebyshev 


MM{tt) 


Fit 


Resolution 


Total 




polynomial order 


binning 


model 






Ml, MeV/c2 


1 O .1/1 

+Z.D /+4 
-0 /-I 


1 n 1 1 1 n 
+0.1 /+0 

-0 / -2 


111 1 1 Q 
+ i.i /+0 

-0 / -0 


+U.3 /+0 

-0 / -0 


+3.0 /+0 
-1.0/ -2 


Ti, MeV 


+1.5 /+9 

-0 /-3 


+ 1.4 /+10 

-0 /-o 


+0.1 /+0 

-0 /-3 


+0 /+0 
-1.2/-1 


+2.1 /+13 
-1.2/-4 


Ms, MeV/c2 


+0 /+1 
-0.4/ -2 


+0 /+2 
-0.6/ -0 


+0 /+0 
-1.4/ -0 


+0 /+0 
-0.3/ -0 


+1.0 /+2 
-1.9/ -2 


Fa, MeV 


+0 /+2 
-5.5/ -2 


+2.1 /+4 

-0 / -0 


+0 1+7 
-0.2/ -0 


+0 /+0 
-0.8/ -1 


+2.1 /+8 
-5.7/ -2 


a 


+0.1 /+0.2 
-0.5/ -0.7 


+0 /+0 
-0/ -0.1 


+0 /+0.4 
-0.1/ -0 


+0 /+0 
-0/ -0 


+0.1 /+0.4 
-0.5/ -0.7 


(p, degree 


+4 /+11 
-4/ -167 


+0 1+5 

-2/ -0 


+0 /+0 
-8/ -74 


+0 /+0 
-l/-0 


+4 /+12 
-9/ -193 



and exclusive yields. The systematic uncertainty in the inclusive yield estimated in Ref. [2| 
is taken into account. When fitting the inclusive MM(7r+7r^) spectra in bins of MM{tt) we 
allow the yield of the reflection to float within the uncertainty of the above measurement. 
Combinatorial background is parameterized by a Chebyshev polynomial function. The order 
of the polynomial varies from 6 to 8 for different MM{tt) bins. The fit results for the hb{2P) 
yield is presented in Fig. [71 

To fit the MM {it) distribution of the hh{2P) signal we use the same fit function as for 
the hb{lP). The result of the fit is shown in Fig. [71 As in the case of the hb{lP), the 
non-resonant amplitude is found to be consistent with zero, 0.11 ± 0.13. Numerical values 
for other parameters are given in Table I VII The confidence level of the fit is 42%. The 
two-resonance hypothesis is favored over that for one resonance (no resonances) at the 2.7 a 
(6.3 0") level. 

The systematic uncertainty is estimated in the same way as for the hf,{lP), the results 
are given in Table IIVI The minimal level at which the two-resonance hypothesis is favored 
over the one-resonance (no-resonance) one for all considered variations is 1.9 cr (4.7 o"). 
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ANGULAR ANALYSIS 



We perform angular analyses to check consistency of the = assignment for the 
^^(lOGlO) and Zf,{10Q50) states and to attempt to discriminate against other hypotheses; 
we consider = 1^, 2+ and 2^. The 0+ (0^) assignment is forbidden by parity conservation 
in Zh —7- T{nS)7T {Z^ — )■ hb{mP)TT) decays. We use the polar angles of the pions (denoted as 
6i for a pion from the Y{5S) decay and 02 for a pion from the Z^ decay), the spatial angle, 
6'7r7r, between the two pions and the angle, 0p, between the plane defined by the pion from 
1^(55*) decay and the beam direction and the plane defined by the two pions. Since the Zj, 
velocity is very small, /3 < 0.02, we neglect its recoil motion and measure all pion momenta 
in the cm. frame. We assume that only the lowest partial wave contributes to the decay in 
the cases where more than one partial wave is possible [3]. 



T(2S')7r+7r- and T(35)7r+7r- final states 

We apply an additional requirement on the angle between the two pion momenta mea- 
sured in the laboratory frame, cos^iab < 0.95, to suppress the background from converted 
photons. 

We first study the angular distribution for the non-resonant component for the 
T(2S')7r^7r^ candidates. The non-resonant region is defined as 10.5GeV/c^ < M(T7r) < 
10.58 GeV/c^. The cos^t^tt values is tightly correlated with the M(7r+7r^), therefore the 
cos^^^TT distribution reflects the details of the M{'K~^7r~) distribution discussed in the Dalitz 
analysis section of the paper. The cos^i, cos ^2 and 0p distributions (after the MM{7i^tt~) 
sideband subtraction and efficiency correction) are all consistent with uniform distributions. 
This indicates that the non-resonant contribution is dominated by an S-wave component. 

At the next step we select the Z;,(10610) [^6(10650)] region 10.602 < M(T(2S)7r) < 
10.638 GeV/c^ [10.648 < M{T{2S)tt) < lO.TGeV/c^] for the T(25)7r+7r- final state and 
10.607 < M(T(3S)7r) < 10.627 GeV/c^ [10.649 < M{T{2S)tt) < 10.676 GeV/c^] for 
the T(3S')7r+7r^ final state. For the T(2S')7r+7r~ candidates we require 0.15GeV^/c'^ < 
M^(7r+7r^) < 0.4GeV^/c^, which considerably suppresses the non-resonant contribution. 
The cos 6*1 , cos 62 and 0p distributions corrected for the reconstruction efficiency are shown 
in Figs E] and El 
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FIG. 8: Angular distributions for the Z;,(10610) signal region. Points with error bars represent 
the yield of the T(25')7r candidates (first row) and T(35')7r candidates (second row) as a function 
of the cos 9i (left column) , cos 02 (middle column) and (j)p (right column) . The open histograms 
represent the fit results for different hypotheses: 1^ (black), 1~ (red), 2+ (green) and 2^ (blue). 
The hatched (cross-hatched) histogram shows the contribution of the non-resonant (combinatorial) 
component. 

We perform a binned maximum likelihood fit to these distributions. The fit function 
is a sum of three components: (a) the Zb signal, (b) the contribution due to the non- 
resonant component and (c) combinatorial background described by a properly normalized 
MM(7r+7r^) sidebands. The terms due to the interference of the Zt signal and 5- wave 
non- resonant part are constant in cos6'i, cos 6*2 and (pp. Therefore, component (b) of the fit 
function is a constant. We estimate its normalization from the linear extrapolation in the 
M(T(2S')7r) for the T(2S')7r+7r^ final state and we assume that it is zero for the T(3S')7r+7r~ 
final state. The non-resonant contribution varies with cos^^,r7r; therefore this variable is not 
used. The components (a) and (b) of the fit function are corrected for efficiency. The only 



fioating parameter in the fit is the normalization of the Zh signal component 



The fit results for various Zb spin-parity assignments are shown in Figs. [HI and [9l We use 
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FIG. 9: Angular distributions for the Zb(10650) signal region. The histogram contents are similar 
to those shown in Fig. [8l 

cos ^2 to discriminate 1^ from 1~ and 1+ from 2^, and we use cos6'i to discriminate 1^ from 
2"*". The probabilities at which the 1~, 2+ and 2^ hypotheses are disfavored compared to 
the l"*" hypothesis are calculated as -\/A21ogL and given in Table |Vl 



TABLE V: The probabilities at which different hypotheses are disfavored compared to the 1"*" 



hypothesis. The data are consistent with the l"*" hypothesis. 





Z5(10610) 


Zfe(10650) 


T(25)7r+^- 


T(35)^+7r- 


/ife(lP)7r+7r- 


T(25)7r+^- 


T(35)7r+7r- 


/ib(lP)7r+7r- 


1" 


3.6(T 


0.3(7 


0.3(7 


3.7(7 


2.6(7 


2.7(7 


2+ 


4.3cj 


3.5cr 


4.3(7 


4.4(7 


2.7(7 


2.1(7 


2" 


2.7a 


2.8(7 


2.9(7 


2.6(7 
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FIG. 10: Points with error bars represent the hh{lP) yield as a function of the cos^i (left column), 
COS02 (middle column) and cos^tttt (right column) for the ^^(10610) (top row) and Z;,(10650) 
(bottom row) regions. Histograms represent the fit results for different hypotheses: 1"*" (black), 
1^ (red), 2+ (green) and 2~ (blue). 

/ib(lP)7r+7r^ final state 

We define the ^^(10610) and Zfe(10650) signal regions as lO.eGeV/c^ < MM{7i) < 
10.63 GeV/c^ and 10.63 GeV/c^ < MM{'k) < 10.67 GeV/c^, respectively. For each region 
we fit the MM(7r+7r^) spectra in bins of cos^^i, cos^^2 and cos^^t^tt to determine the h^^lP) 
signal yield. The /if, (IP) yield corrected for the reconstruction efficiency as a function of 
angular variables is shown in Fig [TOl The results of the fits to the different spin-parity 
hypotheses are superimposed. We use cos^tttt to discriminate l"*" from 2+ and 1"*" from 
2~, and we use cos6'2 to discriminate 1"*" from 1~. The level at which the 1~, and 2~ 
hypotheses are disfavored compared to the 1"*" hypothesis is calculated as a square root of 
the difference of the corresponding values and is shown in Table El 

The values quoted in Table |V] are preliminary. Our procedure to deal with the non- 
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TABLE VI: Comparison of results on ^^(lOGlO) and Zf,(10650) parameters obtained from T(55) — )■ 
T(n5')7r+7r^ (n = 1,2,3) and T(55) — )• hfj(mP)TT^ ir^ (m = 1,2) analyses. Quoted values are in 
MeV/c^ for masses, in MeV for widths and in degrees for the relative phase. Relative amplitude 
is defined as az6(i0650)/aZbi06io- 

Final state T(lS')7r+7r~ T(25)7r+7r- T(35)7r+7r~ /ib(lP)7r+7r- /ib(2P)7r+7r- 

M(Zfe(10610)) 10609 ±3±2 10616 ± 2l^ 10608 ± 2+^ 10605.1 ± 2.2 1^;" 10596 ± 7 1^ 
r(Zfe(10610)) 22.9 ±7.3 ±2 21.1 ±4+^ 12.2 ±1.7 ±4 11.4 +|^+^;^ 16l}[;+f 
M(Zb(10650)) 10660 ± 6 ± 2 10653 ± 2 ± 2 10652 ± 2 ± 2 10654.5 ± 2.5 t\i 10651 ± 4 ± 2 
r(Z6(10650)) 12 ±10 ±3 16.4 ±3.6^^ 10.9 ± 2.6j:^ 20.9l:|-jt|^ ^"^-l^-l 
Rel. amplitude 0.59 ± 0.19l[!;[J| 0.91 ± 0.1lt[!:[J^ 0.73 ± 0.10t[!:J^ 1.8 1 J;? toi 1.3 
Rel. phase, 53 ± QI^Iq -20 ± 18^^^ 6 ± 24^23 i88 +4 255 



resonant contribution is approximate, and we do not take into account mutual cross-feed 
of the ^5(10610) and Z;,(10650) states. However, the angular analyses indicate that the 
l"*" hypothesis for both Z;,(10610) and Z;,(10650) provides the best description of angular 
distributions compared to all other hypotheses with J < 2. 



DISCUSSION AND CONCLUSIONS 



In conclusion, we have observed two charged bottomonium-like resonances, the ^^(10610) 
and ^{,(10650), with signals in five different decay channels, T{nS)TT^ {n = 1,2,3) and 
hblmP)^^ (m = 1,2). Parameters of the resonances as measured in different channels are 
summarized in Table IVIi All channels yield consistent results as can be seen in Fig. [HI A 
simple weighted averages over all five channels give M[Z;,(10610)] = 10608.4 ± 2.0MeV/c^, 
F[Z;,(10610)] = 15.6 ±2.5 MeV and M[Z;,(10650)] = 10653.2 ± 1.5 MeV/c^, F[Z5(10650)] = 
14.4 ± 3.2 MeV, where statistical and systematic errors are added in quadrature. 

Charged bottomonium-like resonances cannot be simple bb combinations. The measured 
masses of these new states exceed by only a few MeV/c^ the thresholds for the open beauty 
channels B*B (10604.6 MeV) and B*B* (10650.2 MeV). This "coincidence" is suggestive of 
"molecular" states whose their structure is determined by the strong interaction dynamics 
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FIG. 11: Comparison of Z;,(10610) and Z;,(10650) parameters obtained from different decay chan- 
nels. The vertical dotted lines indicate B*B and B*B* thresholds. 

of B* B and B* B meson pairs. 

The widths of both states are similar and are of the order of 15MeV/c^. The Zfe(10610) 
production rate is similar to the ^^(10650) production rate for every decay channel. Their 
relative phase is consistent with zero for the final states with the T{nS) and consistent with 
180 degrees for the final states with hi,{mP). 

The T(5S') — )■ hb{mP)TT~^7!-~ decays seem to be saturated by the Zfe(10610) and ^^(10650) 
intermediate states; this decay mechanism is responsible for the high rate of the T{5S) — ?• 
hb{mP)TT^TT~ process measured recently by the Belle Collaboration. 

Analysis of angular distributions for charged pions favors the = l"*" spin-parity assign- 
ment for both ^6(10610) and Z;,(10650). Since the T{5S) has negative G-parity, Zf, states 
will have opposite G-parity due to emission of the pion. 
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